! Bessel Y_1(x) function in double precision
!
      function dbesy1(x)
      implicit real*8 (a - h, o - z)
      dimension a(0 : 15), b(0 : 111), c(0 : 125), d(0 : 51)
      parameter (pi4 = 0.78539816339744830962d0)
      data (a(i), i = 0, 15) / 
     &    -0.00000000000009464990d0, 0.00000000002141432796d0, 
     &    -0.00000000359779443470d0, 0.00000043173541397185d0, 
     &    -0.00003453883313621946d0, 0.00165786399054057257d0, 
     &    -0.03978873577297383381d0, 0.31830988618379067154d0, 
     &    -0.00000000005571931463d0, 0.00000000893098810365d0, 
     &    -0.00000099267351748541d0, 0.00007164268731546281d0, 
     &    -0.00295530533604595764d0, 0.05434868816050718549d0, 
     &    -0.19605709064623884320d0, -0.63661977236758134367d0 / 
      data (b(i), i = 0, 13) / 
     &    0.00000000002553685555d0, 0.00000000007384080449d0, 
     &    -0.00000000029605446927d0, -0.00000000520817391060d0, 
     &    -0.00000003676554001376d0, -0.00000015326525336270d0, 
     &    -0.00000013967947147386d0, 0.00000388730181348638d0, 
     &    0.00004052240780599446d0, 0.00025133005806772436d0, 
     &    0.00109380629039004777d0, 0.00298157276348996763d0, 
     &    0.00047188983238582062d0, -0.19643830264444486066d0 / 
      data (b(i), i = 14, 27) / 
     &    0.00000000001056916126d0, 0.00000000062160122556d0, 
     &    0.00000000385741197273d0, 0.00000000749146940600d0, 
     &    -0.00000006306874425799d0, -0.00000068698336577990d0, 
     &    -0.00000329785491054157d0, -0.00000557893716372629d0, 
     &    0.00004634569030169549d0, 0.00049086790937752985d0, 
     &    0.00256488729380268318d0, 0.00822571374211646414d0, 
     &    0.01094512145046992425d0, -0.19159562910398869317d0 / 
      data (b(i), i = 28, 41) / 
     &    -0.00000000040211480479d0, -0.00000000229721816531d0, 
     &    0.00000000062931855684d0, 0.00000006509879184236d0, 
     &    0.00000031589903600713d0, -0.00000003093448289001d0, 
     &    -0.00000845532004628483d0, -0.00004930719004112743d0, 
     &    -0.00009865914153735679d0, 0.00047059803263198826d0, 
     &    0.00472271367216891688d0, 0.01915518502418162917d0, 
     &    0.03722391361384158740d0, -0.16933230788966878110d0 / 
      data (b(i), i = 42, 55) / 
     &    0.00000000110728466147d0, 0.00000000205608441891d0, 
     &    -0.00000002536199750468d0, -0.00000011751665325851d0, 
     &    0.00000034045343850770d0, 0.00000436083390507752d0, 
     &    0.00000834725302600388d0, -0.00007154971448642437d0, 
     &    -0.00051956165332380432d0, -0.00100811886858456098d0, 
     &    0.00436854257164588225d0, 0.03425644258611852849d0, 
     &    0.09073971172014934990d0, -0.10791597118914788657d0 / 
      data (b(i), i = 56, 69) / 
     &    -0.00000000117749768422d0, 0.00000000168403859797d0, 
     &    0.00000003471030171126d0, -0.00000006227816248757d0, 
     &    -0.00000118405395385178d0, 0.00000014478514453927d0, 
     &    0.00003630744090017806d0, 0.00010481065165516454d0, 
     &    -0.00052483993430622472d0, -0.00407014291708523893d0, 
     &    -0.00574311514714191566d0, 0.03534768318503458091d0, 
     &    0.16539360869993020163d0, 0.01986240425870529979d0 / 
      data (b(i), i = 70, 83) / 
     &    0.00000000040080505520d0, -0.00000000336044652853d0, 
     &    -0.00000000092387130705d0, 0.00000017840823064933d0, 
     &    -0.00000022643988836534d0, -0.00000808896250262948d0, 
     &    -0.00000169403516892975d0, 0.00026617475073041589d0, 
     &    0.00072562705080786576d0, -0.00399208199997748377d0, 
     &    -0.02399816920965700841d0, -0.00925574744557542440d0, 
     &    0.20082814022820666990d0, 0.21040564991987672576d0 / 
      data (b(i), i = 84, 97) / 
     &    0.00000000006692975871d0, 0.00000000016722027050d0, 
     &    -0.00000001433097981251d0, 0.00000005455077276306d0, 
     &    0.00000104954046999028d0, -0.00000339572219821681d0, 
     &    -0.00005541577161092338d0, 0.00004344113487321887d0, 
     &    0.00184798594297387512d0, 0.00301041508468150450d0, 
     &    -0.02789739141907399514d0, -0.09421628519129373841d0, 
     &    0.09950212445529771576d0, 0.37496975838746404080d0 / 
      data (b(i), i = 98, 111) / 
     &    0.00000000004528583204d0, 0.00000000066861092459d0, 
     &    -0.00000000924648168570d0, -0.00000008002058622489d0, 
     &    0.00000087578676979688d0, 0.00000628456150216647d0, 
     &    -0.00004276239654834506d0, -0.00034630842915911896d0, 
     &    0.00089431202169021474d0, 0.01086380107342193819d0, 
     &    0.00145563493724746238d0, -0.14193383817937981478d0, 
     &    -0.15148562296986494725d0, 0.35720232689066512346d0 / 
      data (c(i), i = 0, 13) / 
     &    0.00000000032016204980d0, -0.00000000177431944127d0, 
     &    0.00000001367520951871d0, -0.00000006171755212447d0, 
     &    -0.00000043602779932091d0, 0.00000179827173442894d0, 
     &    0.00004990624965075363d0, -0.00018618453816525107d0, 
     &    -0.00206578613620093339d0, 0.00640292258731756063d0, 
     &    0.04663627606907821596d0, -0.11400070635033154332d0, 
     &    -0.26159330264498333979d0, 0.30099732306965462346d0 / 
      data (c(i), i = 14, 27) / 
     &    -0.00000000000059774708d0, -0.00000000038827648808d0, 
     &    0.00000000527269271362d0, 0.00000002971059814012d0, 
     &    -0.00000052701406514027d0, -0.00000320369808906187d0, 
     &    0.00004475144442564492d0, 0.00016907703884518866d0, 
     &    -0.00209937055264007564d0, -0.00491023437935628693d0, 
     &    0.04967175540701545140d0, 0.04195559522923170315d0, 
     &    -0.33516091307165838036d0, -0.02375823895638961835d0 / 
      data (c(i), i = 28, 41) / 
     &    0.00000000000283380103d0, -0.00000000044715450235d0, 
     &    0.00000000017695199610d0, 0.00000006053379840070d0, 
     &    -0.00000002888550814408d0, -0.00000593780185020349d0, 
     &    0.00000514027455732102d0, 0.00035667482807472357d0, 
     &    -0.00038132729962607730d0, -0.01158750857180657934d0, 
     &    0.01376562658346365929d0, 0.14388461069626297876d0, 
     &    -0.13107454661755534290d0, -0.27409127395927545296d0 / 
      data (c(i), i = 42, 55) / 
     &    0.00000000002683697058d0, -0.00000000024302161327d0, 
     &    -0.00000000430151236236d0, 0.00000003557058903199d0, 
     &    0.00000049328952592926d0, -0.00000365732389065382d0, 
     &    -0.00003642356546467712d0, 0.00023636175253326716d0, 
     &    0.00154549935517002108d0, -0.00837084637397210913d0, 
     &    -0.02941202173717025197d0, 0.11713664042453548043d0, 
     &    0.15186375421302413107d0, -0.25912851048611625179d0 / 
      data (c(i), i = 56, 69) / 
     &    0.00000000002914434604d0, 0.00000000015185121893d0, 
     &    -0.00000000488084673389d0, -0.00000001936460341189d0, 
     &    0.00000057968735781303d0, 0.00000158978381365925d0, 
     &    -0.00004521858314881102d0, -0.00007426575562367423d0, 
     &    0.00206293925392279862d0, 0.00143934355204571641d0, 
     &    -0.04414786410004317489d0, -0.00317227451807568106d0, 
     &    0.27328377353032129600d0, -0.02616867939853747003d0 / 
      data (c(i), i = 70, 83) / 
     &    0.00000000000644550655d0, 0.00000000040249587279d0, 
     &    -0.00000000123789833341d0, -0.00000005582997578770d0, 
     &    0.00000016975984093144d0, 0.00000524003568623094d0, 
     &    -0.00001539871212028104d0, -0.00030372172637669852d0, 
     &    0.00082284448843625174d0, 0.00923614667011392611d0, 
     &    -0.02069442879928366384d0, -0.10834973445898588438d0, 
     &    0.14982326837249145873d0, 0.20317989938720766823d0 / 
      data (c(i), i = 84, 97) / 
     &    -0.00000000002115691215d0, 0.00000000029752035314d0, 
     &    0.00000000334739862269d0, -0.00000004304625868398d0, 
     &    -0.00000036733494871207d0, 0.00000425328270524107d0, 
     &    0.00002584872967637478d0, -0.00026244994622929809d0, 
     &    -0.00102234403086242668d0, 0.00863185023199066935d0, 
     &    0.01816610124346617164d0, -0.11151661591132254182d0, 
     &    -0.08978791805571149950d0, 0.23370422835726857838d0 / 
      data (c(i), i = 98, 111) / 
     &    -0.00000000002957977086d0, -0.00000000006150420279d0, 
     &    0.00000000488089472193d0, 0.00000000624218182859d0, 
     &    -0.00000056561995267846d0, -0.00000032077350100271d0, 
     &    0.00004278569848964830d0, -0.00000053069698280036d0, 
     &    -0.00187141538785600230d0, 0.00073172821883575326d0, 
     &    0.03832854206202077374d0, -0.01874044416229564036d0, 
     &    -0.23027059405144880640d0, 0.05794254714300082167d0 / 
      data (c(i), i = 112, 125) / 
     &    -0.00000000001203179184d0, -0.00000000035509887475d0, 
     &    0.00000000213620090322d0, 0.00000004813181550803d0, 
     &    -0.00000026815627351126d0, -0.00000439191902440609d0, 
     &    0.00002211921535574542d0, 0.00024593505480524700d0, 
     &    -0.00106158930833741898d0, -0.00722666167492552731d0, 
     &    0.02396424370522244155d0, 0.08278313570069734843d0, 
     &    -0.15890724632166919294d0, -0.15383825653750118007d0 / 
      data (d(i), i = 0, 12) / 
     &    -1.272346002224188092d-14, 3.370464692346669075d-13, 
     &    -1.144940314335484869d-11, 6.863141561083429745d-10, 
     &    -9.491933932960924159d-8, 5.301676561445687562d-5, 
     &    0.1628675039676399740d0, -3.652982212914147794d-13, 
     &    1.151126750560028914d-11, -5.165585095674343486d-10, 
     &    4.657991250060549892d-8, -1.186794704692706504d-5, 
     &    1.562499999999994026d-2 / 
      data (d(i), i = 13, 25) / 
     &    -8.713069680903981555d-15, 3.140780373478474935d-13, 
     &    -1.139089186076256597d-11, 6.862299023338785566d-10, 
     &    -9.491926788274594674d-8, 5.301676558106268323d-5, 
     &    0.1628675039676466220d0, -2.792555727162752006d-13, 
     &    1.108650207651756807d-11, -5.156745588549830981d-10, 
     &    4.657894859077370979d-8, -1.186794650130550256d-5, 
     &    1.562499999987299901d-2 / 
      data (d(i), i = 26, 38) / 
     &    -6.304859171204770696d-15, 2.857249044208791652d-13, 
     &    -1.124956921556753188d-11, 6.858482894906716661d-10, 
     &    -9.491867953516898460d-8, 5.301676509057781574d-5, 
     &    0.1628675039678191167d0, -2.185193490132496053d-13, 
     &    1.048820673697426074d-11, -5.132819367467680132d-10, 
     &    4.657409437372994220d-8, -1.186794150862988921d-5, 
     &    1.562499999779270706d-2 / 
      data (d(i), i = 39, 51) / 
     &    -4.740417209792009850d-15, 2.578715253644144182d-13, 
     &    -1.104148898414138857d-11, 6.850134201626289183d-10, 
     &    -9.491678234174919640d-8, 5.301676277588728159d-5, 
     &    0.1628675039690033136d0, -1.755122057493842290d-13, 
     &    9.848723331445182397d-12, -5.094535425482245697d-10, 
     &    4.656255982268609304d-8, -1.186792402114394891d-5, 
     &    1.562499998712198636d-2 / 
      if (x .lt. 0.85d0) then
          t = x * x
          y = (((((((a(0) * t + a(1)) * t + 
     &        a(2)) * t + a(3)) * t + a(4)) * t + 
     &        a(5)) * t + a(6)) * t + a(7)) * x
          y = (((((((a(8) * t + a(9)) * t + 
     &        a(10)) * t + a(11)) * t + a(12)) * t + 
     &        a(13)) * t + a(14)) * t + a(15)) / x + 
     &        y * log(x)
      else if (x .lt. 4.15d0) then
          v = 4 / x
          t = x - v
          k = int(t + 4)
          t = t - (k - 3.5d0)
          k = k * 14
          y = (((((((((((((b(k) * t + b(k + 1)) * t + 
     &        b(k + 2)) * t + b(k + 3)) * t + b(k + 4)) * t + 
     &        b(k + 5)) * t + b(k + 6)) * t + b(k + 7)) * t + 
     &        b(k + 8)) * t + b(k + 9)) * t + b(k + 10)) * t + 
     &        b(k + 11)) * t + b(k + 12)) * t + b(k + 13)) * v
      else if (x .lt. 12.5d0) then
          k = int(x)
          t = x - (k + 0.5d0)
          k = 14 * (k - 4)
          y = ((((((((((((c(k) * t + c(k + 1)) * t + 
     &        c(k + 2)) * t + c(k + 3)) * t + c(k + 4)) * t + 
     &        c(k + 5)) * t + c(k + 6)) * t + c(k + 7)) * t + 
     &        c(k + 8)) * t + c(k + 9)) * t + c(k + 10)) * t + 
     &        c(k + 11)) * t + c(k + 12)) * t + c(k + 13)
      else
          v = 24 / x
          t = v * v
          k = 13 * (int(t))
          y = ((((((d(k) * t + d(k + 1)) * t + 
     &        d(k + 2)) * t + d(k + 3)) * t + d(k + 4)) * t + 
     &        d(k + 5)) * t + d(k + 6)) * sqrt(v)
          theta = (((((d(k + 7) * t + d(k + 8)) * t + 
     &        d(k + 9)) * t + d(k + 10)) * t + d(k + 11)) * t + 
     &        d(k + 12)) * v - pi4
          y = y * (-cos(x + theta))
      end if
      dbesy1 = y
      end
!
